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We study the probability distribution of residence time of a grain at a site, and its total residence 
time inside a pile, in different ricepile models. The tails of these distributions are dominated by the 
grains that get deeply buried in the pile. We show that, for a pile of size L, the probabilities that 
the residence time at a site or the total residence time is greater than t, both decay as l/t(\nt) x for 
L u Ct< exp(L 7 ) where 7 is an exponent > 1, and values of x and ui in the two cases are different. 
In the Oslo ricepile model we find that the probability that the residence time Ti at a site i being 
greater than or equal to t, is a non-monotonic function of L for a fixed t and does not obey simple 
scaling. For model in d dimensions, we show that the probability of minimum slope configuration 
in the steady state, for large L, varies as exp( — KL d+2 ) where k is a constant, and hence 7 = d + 2. 



1. INTRODUCTION 

Granular materials have drawn a lot of attention due 
to their complex flow behaviour under different driving 
conditions pjj. Slowly driven pile of sand grains serve as 
a prototype for self-organized criticality (SOC) 0. Al- 
though SOC was not seen in experiments on piles of sand 
0, but experiments on piles of long grained rice have 
shown evidence of power law distribution of avalanche 
sizes 0,0,0- Historically, perhaps because of its relation 
to earthquake phenomena, studies of sandpiles have 
generally focused on the distribution of avalanche sizes. 
There are only a few theoretical studies of other interest- 
ing quantities such as the distribution of total residence 
times of grains in piles, even though the experimental 
studies by the Oslo group 0, using coloured tracer 
grains are now almost a decade old. 

In this paper we consider the probability distribution 
of residence times of grains at a site, and of their 
total residence times in the pile, in critical slope type 
slowly-driven sandpile (equivalently ricepile) models. 
In these models, these residence time distributions 
are qualitatively different from the critical height type 
models. In critical height models, the distribution 
decays exponentially with average total residence time 
equal to average active mass in the pile 0- In critical 
slope models, there is a possibility that the grain gets 
buried very deep in the pile, and then takes a long 
time to come out. We shall show that this makes the 
cumulative probability distribution of these residence 
times to have a characteristic 1/t decay for large times 
t, (modified by a logarithmic multiplicative correction 
factor), independent of details of the toppling rules, and 
of the dimensionality of the system. 

There have been some numerical and analytical studies 
of these distributions earlier. Frette has proposed a the- 
oretical model, called the Oslo ricepile model which 
seems to reproduce the phenomenology of the ricepile ex- 
periment well. In the experimental studies of ricepiles, 
Christensen et. al. 0, 0] estimated that average total 



residence time of grains in a pile of size L varies as L v , 
where they estimated v — 1.5 ±0.2. From numerical sim- 
ulations of the Oslo model for systems of size L < 1600, 
the exponent characterizing the power law decay of the 
probability density of total residence times at large times 
was estimated as 2.2 ± 0.1. Boguna and Corral jlOj , and 
Carreras et. al. have used a continuous-time ran- 
dom walk model of the motion of grains, with long trap- 
ping times and a power-law distribution of step sizes, 
to explain the anomalous diffusion of tracer grains. In 
earlier papers 0, > we have studied the total residence 
time distribution in the critical height type sandpile mod- 
els with both deterministic and stochastic toppling rules. 
We reduced the problem to a diffusion problem of a single 
particle in a medium with space dependent jump rates 
and showed that the distribution of the total residence 
time does not have any power law tail. We also obtained 
the non-universal scaling form of the distribution which 
depends on the probability distribution of where grains 
are added into the pile. 

In this paper we study the distribution of total resi- 
dence times of grains in the pile, and also of residence 
time at a site, in a class of critical slope type sandpile 
models. We define the residence time Ti at a site i is the 
time spent by a grain at the this site, measured in units 
of the time interval between successive addition of grains. 
The total residence time T is defined similarly. We show 
that the probability of the residence time at a site or the 
total residence time in the pile, being greater than or 
equal to t, decays as l/t(\nt) s for a very wide range of 
t. The upper cutoff in both the distributions scales with 
system size L as exp(KL 7 ) where 7 is an exponent > 1 
and k is a positive constant. 

For the Oslo ricepile model, we find an unexpected be- 
haviour in the cumulative probability that a grain staying 
at a site i at least upto time t, is not a monotonically in- 
creasing function of system size L. We argue that this 
implies the cumulative probability distribution function 
ProbL{T\ > t) cannot have a simple finite size scaling 
form. We show that 7 = d + 2 in d dimensions for this 
model. 

Plan of the paper is as follows. In section 2 we define 
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the four models studied in this paper. In section 3 we 
present the simulation results for the residence time T\ 
at site 1 for the Id Oslo ricepile model and explain the 
non-monotonic behaviour of the cumulative distribution 
ProbL (Ti > t) with L by relating the residence times of 
grains at the site 1 to the statistical properties of height 
fluctuation at that site. We also explain the origin of 
multiplicative logarithmic correction factor appearing in 
the 1/t decay of ProbL{Ti > t). In section 4 we discuss 
the 1/t power law form of ProbhiT > t), where T is the 
total residence times, for large t in the Id Oslo model 
and show that this also has a multiplicative logarithmic 
correction. In section 5 we argue that the probability 
of minimum slope configuration occurring in the steady 
state of the Id Oslo ricepile model, scales with system 
size L as exp(— kL 3 ) where k is some positive constant. 
In section 6 we present our simulation results for other 
models and show that in all cases the cumulative dis- 
tributions is qualitatively similar to the Id Oslo ricepile 
model. The last section contains a summary and some 
concluding remarks. 

2. DEFINITION OF THE MODELS 

We consider general critical slope type sandpile models 
where the configurations are specified by integer height 
variables h(x), i.e., number of grains, at any site a? of a 
finite d-dimensional lattice. Whenever height difference 
between two adjacent sites is greater than a threshold 
value, some specified number of grains are transferred to 
the neighbouring sites. Piles are driven by adding grains, 
one at a time, at a fixed, or at a randomly chosen site. 
Grains are added only when there are no unstable sites 
left in the system, and can leave the pile from the bound- 
ary. We update all unstable sites in parallel. We have 
studied four different models both in one and two dimen- 
sions : the Oslo ricepile model and it's 2d generalization, 
local limited model and it's variation. We now define the 
precise rules of these four models. 



(rightmost end), one grain goes out of the system. Grains 
are added only at site 1. 

The Id Oslo ricepile model has an abelian property 
The final height configuration does not depend on 
the order we topple the unstable sites. After addition of 
total L(L + 1) grains, the pile reaches the critical steady 
state [l3J. Since we have chosen the values of z c to be 
1 or 2, in the steady state height profile fluctuate be- 
tween slope 1 and 2. For number of sites L, the number 
of possible configurations in the critical states are expo- 
nentially large, approximately ( 3+ 2 ^ f° r large L 
[l^ . The probabilities of various configurations in the 
steady state differ from one another by many orders of 
magnitude unlike the BTW model 

Model-B : 2d generalization of the Oslo model. 

The Oslo model defined above can easily be general- 
ized to two dimensions. We take a triangular region of 
a square lattice, the sites of which are indexed by 
with i,j > 1 and i + j < L + 1. The height of the pile 
at site is denoted by h(i,j). Whenever the height 
difference between site and any of it's neighbouring 
sites exceeds a critical value z c (i,j), assigned to the site 
there is a toppling at site and one grain is 

transferred from this site to the lower neighbouring site 
towards the unstable direction. If there are more than 
one unstable directions, grain is transferred towards the 
greatest slope. If the two directions have equal slope val- 
ues, one grain is transferred randomly towards any one 
of these two directions. Whenever there is a toppling at 
a site (i, j), z c (i,j) is reset randomly, independent of the 
history, to either 2 or 1 with probability p or q respec- 
tively, where p + q = 1. One grain is lost, whenever there 
is a toppling at the boundary sites i.e., along i+j = L + l 
line. The model defined above in two dimensions is not 
abelian because final stable configuration depends on the 
order we topple the unstable sites. Grains are added only 
at the corner site (1,1). 



Model-A: The Oslo ricepile model. 

The Oslo ricepile model || is defined as follows. We 
consider a one dimensional ricepile, which is specified 
by an integer height variable hi at each site i of a one- 
dimensional lattice, with 1 < i < L. The slope at site i 
is defined to be hi — hi + \. Whenever the slope z\ at any 
site i is higher than a critical value z c j, the site becomes 
unstable and one grain from the unstable site goes to the 
right neighbour, i.e., hi — > hi — 1 and — > + 1. 
Whenever there is a toppling at site i, z c% i is randomly, 
independent of the history, reset to one of the two val- 
ues, 1 and 2, with probability q and p respectively, where 
p + q = 1 . Whenever there is a toppling at site i = L 



Model-C : The local limited model. 

The local limited model [l5j is a one dimensional model 
defined as follows. The slope Zi is defined as i.e. , Zi — 
hi — hi + \. Whenever value of the slope Zi at any site 
i is higher than a critical value z c , which we choose to 
be 2, the site becomes unstable and two grains from the 
unstable site goes to the right neighbour, i.e., hi — > hi — 2 
and hi + i — ► h i+ i + 2. Slope at any site may be negative 
in the local limited model. Whenever there is a toppling 
at site i = L (rightmost end), two grains go out of the 
system simultaneously. 

Grains are added uniformly everywhere. This model 
is also not abelian. It is easy to see that in this case 
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FIG. 1: Rice pile of size L = 5 after addition of 100 grains. All 
grains are numbered whenever added in the pile. Minimum 
slope is denoted by the thick line. 



the maximum and the minimum slopes are 2 and 1 re- 
spectively. Total number of recurrent configurations in 
the steady state can be determined exactly, and varies as 
-p72 for large L ^(| . It is known that the probabilities of 
occurrence of various configurations in the steady state 
is not equal, and may differ from one another by many 
orders of magnitude. 



other (Fig. 1). Whenever a grains is added at a site, it 
sits on the top of the stack. When one unstable grain 
leaves the stack, it is taken from the top of the stack. In 
model-C, when two grains leave a site, we first take out 
the topmost grain from the site and put it on the top of 
right nearest neighbour stack, then we take the second 
unstable grain and put it on the top of the first grain at 
right nearest stack. In model-D, we transfer the unstable 
grain, second from the top, to the right nearest neigh- 
bour and transfer the topmost one to next to the right 
nearest neighbour. 

If a particular grain n enters a site i at time t in (i,n) 
and leaves the site at time t out (i,n), it's residence time 
Tj(n) at site i is defined as the time spent by the grain at 
the site 2, i.e., Tj(n) = t out {i,n) —ti n (i,n). The residence 
time of the n th grain, T(n), is the total time spent by 
the grain inside the pile. For a directed ricepiles in one 
dimension where grains move only in one direction and 
by one step in each toppling, the residence time T(n) 
equals to ^2i=iTi{n) {e.g. in model-A and model-C). 
We define the function ProbhiTj > t) as the probability 
that a new grain added in the steady state of the pile will 
have a residence time at site j is greater than or equal 
to t, and Prob^iT > t) as the probability that its total 
residence time in the pile is greater than or equal to t. 
Clearly, we have Prob L (Tj > 0) = Prob L {T > 0) = 1. 



RESIDENCE TIMES AT THE FIRST SITE IN 
THE OSLO RICEPILE MODEL 



Model-D : Model with non-nearest neighbour 
transfer of grains. 

Model-D is a variation of the model-C |l5|. Whenever 
value of the slope Zi > 2 at any site i, the site becomes 
unstable and two grains from the unstable site are 
transferred to the right, one grain transferred to site 
i + 1 and the other one transferred to site i + 2, i.e., 
hi — > hi - 2, h i+ i — > h l+ i + 1 and h i+2 — ► h i+2 + 1. 
If there is a toppling near the right boundary, grain 
goes out of the pile. The order we relax unstable sites 
matters. The grains are added uniformly everywhere. 
The local slope can be negative as in model-C. The 
minimum and maximum slope in this model are also 1 
and 2 respectively. 

The pile in all four cases is driven slowly, by adding 
one grain per unit time, starting with the initial config- 
uration of height zero at all sites. We assume that the 
time interval between addition of two grains is chosen 
long enough so that all avalanche activity has died be- 
fore a new grain is added. The grain added at time n 
will be labeled by the number n. We think of the grains 
at a particular site as stacked vertically, one above the 



The qualitative behaviour of distributions Prob^Ti > 
t) for i — 1 can be seen in the simulation results shown 
in Fig. 2 and Fig. 3. We have done our simulations for 
p = q = I and different system sizes, L = 20,25,35 and 
50. We averaged the data for a total 10 9 grains added in 
the pile for each L. Fig. 2 shows the plot of ProbL(T\ > 
t) versus time t for different values of L. Interestingly, 
various curves for different L have steps like structures. 
The curves for different values of L cross each other many 
times. The unusual non-monotonic behaviour is not an 
artifact of statistical fluctuations. The statistical errors 
in the data are much smaller than the step sizes except in 
the tail region (i.e, t 10 6 ). The crossing of the curves 
for the cumulative probabilities persists for quite large 
system sizes also. In Fig. 3, we have plotted Prob^iTx > 
t) versus t for two much bigger system sizes, L = 100 and 
L = 200. We see that in this case also the probability 
that a grain remains in the pile of size L = 100, for 
time greater than or equal to 6 x 10 5 , is higher by a 
factor 1.8 than for a pile, two times larger size L = 200. 
Steps like structures are not log periodic as the height 
and width of a step in each curve increases when going 
down the curve even on the log scale. The existence of 
several steps, whose positions and logarithmic widths are 
different for different L's, implies that simple finite size 
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FIG. 2: The cumulative probability ProbhiTi > t) versus 
time t for lattice sizes L = 20, 25, 35 and 50 in the Id Oslo 
ricepile model. Total 10 9 grains were added. 



FIG. 4: Scaling collapse of various probability distributions 
Profoi(A/ii) where A/ii is the deviation of height at site 
1 about it's average value, for different system sizes, L — 
100, 200 and 400 for the Id Oslo ricepile model. 
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FIG. 3: The cumulative probability Probz,(Ti > t) versus 
time t for lattice sizes L = 100 and L = 200 in the Id Oslo 
ricepile model. 



scaling cannot hold in this case. 



Relationship between residence times Xi and height 
fluctuations. 

We can understand the residence time distribution of 
grains at any site in terms of the fluctuation of height at 
that site. The hi(i) be the height of the pile at a site i just 
after the t th grain has been added. This is a stochastic 
process and, in the steady state, it fluctuates in time 
between a upper bound, h max , and a lower bound, h min . 
In case of the height fluctuation at site 1, h max = 2L and 
hmin — L. The height hi(t) at the site 1 has a stationary 
probability distribution which is sharply peaked near its 
average value hi, and has the width a h which is standard 
deviation of the fluctuation of height hi. In the steady 



state, the average value of hi varies as L, and the width 
cr h varies as L Wl , where exponent wi < 1. For large L, 
the probability distribution of hi has a scaling form as 
given below. 

Prob i (/ ll )= L- Ul 9(^j^) (1) 

In Fig. 4 we have shown a scaling collapse of various 
probability distribution of height at site 1, Probi,(Ahi), 
where Ahi = hi — hi, for various values of system sizes, 
L = 100, 200 and 400 in the Id Oslo ricepile model. We 
get a good collapse using the scaled variable Ah/L^ 1 
where u>i ~ 0.25. Here the scaling function g(x) is nearly 
Gaussian for x near zero. But very large deviations of hi 
from the mean value are not well-described in the Gaus- 
sian approximation. Later we shall argue that in the Oslo 

i 

model scaling function g(x) varies as exp(— \x\ ) for 

3 

x » 1 and it varies as exp(— \x\ ) for x <C —1. 

Let us consider variation of height hi at the first site 
with time t shown schematically in Fig. 5. Note that 
hi(t) is piecewise constant line segments, with possible 
jumps at the integer time t. The value of hi(t) at time 
t is denoted by y-coordinate of the line segment which 
is just at the right of the coordinate t, e.g., hi(0) = 33, 
hi(l) — 34, etc. A grain added at time t, when the height 
at the first site is hi(t — 1), leaves the site at time t', we 
must have hi{t') < hi(t), and h x {t") > h x {t), for all t" 
satisfying t < t" < t' . As an example, for the time series 
of hi(t) shown in Fig. 5, the grain added at t — 14 stays 
at site 1 upto time t = 41 and then goes out of the site 1 
at t = 42 (i.e., just after addition of the 42nd grain), and 
so Ti(14) = 28. As hi(U) = M 12 ), the grain added at 
t = 12 comes out immediately, and hence Ti(12) = 0. 
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FIG. 5: Time fluctuation of height at first site is plotted as a 
function of time in the Id Oslo ricepile model for L = 20. The 
horizontal line is at hi — 33 which is the maximum probable 
height. The first and second vertical lines are at t = 14 and 
t — 42 respectively. 



Let ProbL{Ti > t\h{) be the conditional probabil- 
ity that a grain stays at site 1 for time greater than t, 
given that it was added when the height was h%. Since 
ProbL(hi) is the probability that height was hi when the 
grain was added, we have the following, summing over all 
possible values of hi. 



Prob L {Ti >t) 



£ 

hi—hfn- 



Prob L {hi)Prob L {Ti > t\hi) 



(2) 

But Probh{Ti > t\hi) can also be written as the con- 
ditional probability that the height of the pile at site 1 
would remain above hi for an interval > t, given that 
the height is hi in the steady state. This probability 
can be calculated from the general theory of Markov 
chains as the probability of first return to a height less 
than or equal to hi, given that we start with height hi 
in the steady state, and add one grain per unit time. 
The probability that no return has occurred up to time 
t decreases as exp[— X(hi)t] for large t, where X(hi) is 
the largest eigenvalue of the the reduced Markov matrix, 
with rows and columns corresponding to configurations 
with heights at site 1, below or equal to hi, removed 
0,0]. While it is not very easy to calculate X(hi) ex- 
actly, clearly it decreases as hi decreases. For hi = h max , 
it is +oo as the height at the site cannot be higher than 
hmax and Ti must always be zero. Also it is very small 
for hi near h m i n , as the pile returns to very low values 
of hi only rarely. 

For large t, in the sum in r.h.s. of Eq.(J2J), only terms 
with hi near h m i n make a significant contribution. In 
this case, it is a reasonable approximation to replace the 



function Prob^Ti > t\hi) by a simple exponential, with 
\{h{) = {Ti)^- Thus we write, for large t, 



Prob L (Ti > t\hi) ~ exp(-t/<Ti) hl ) 



(3) 



It is easy to write the conditional expectation value 
of the residence time at the first site, {Ti)h 1 , given that 
the grain was added at the height hi in terms of the 
stationary probability distribution Pro6i(/ii) exactly as 



l hi 



Prob ^{height > hi) 
PiProb L (hi) 



(4) 



where pi is the probability of adding a grain at site i = 1 . 
When we add grains only at first site, pi = 1 and when 
we add grains uniformly everywhere, pi = l/L. 

Proof : Define an indicator function rj n t = 1, if the 
n th grain is at height hi at time t, and zero otherwise. 
Clearly, the sum of r\ n ^t over t is the residence time of 
n th grain at height hi- Then, averaging over n we get 
the mean residence time. But the sum of r\ n t over n and 
t both gives a contribution whenever there is a grain at 
height hi, and hence is equal to N ProbL^height > hi) 
where N is total number of grains added and N is very 
large. Dividing this sum by average number of grains 
added at height hi, which is equal to piN Prob , we 
get (Ti) hl . Hence, Eq.(@} follows. 

We substitute this estimate of (Ti)/^ in Eq. 0. We 
note that for large t, the terms in the summation that 
contribute significantly correspond to hi near h m i n . For 
these values of hi, Prob(height > hi) is nearly 1, and 
(Ti)/ ll may be replaced, with small error, by 1/ ' Probi,(hi) 
(see Eq. . Then Eq. can be approximately written 
as given below. 



Prob L {Ti>t)~ 2 Prob L (hi)e- tProH{hl) (5) 

hi — hmin 

Thus, the distribution of residence times Ti can 
be expressed in terms of the probability distribution 
ProbL(hi) of height hi. 



Behaviour of ProbL(T\ > t) for large t. 

Now we shall use the knowledge of the behaviour of 
ProbL(hi) to explain the step-like structures in the dis- 
tribution function ProbL^Ti > t). 

For hi <C hi, the probability distribution of height 
ProbL(hi) falls very rapidly. Actually, it will be ar- 
gued in section 5 that for hi hi the ratio Probhihi — 
\)jProbh{hi) is of order exp(— aL 2 ) where a is a con- 
stant and hence is very much less than 1. The values 
of ProbLihi) for different his could differ by several or- 
ders of magnitude from each other, if hi is sufficiently 
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FIG. 6: The cumulative probability ProbhiTi > t) versus 
time t for lattice size L = 20. 



near h min . Now in the interval of 1/ 'Prob^hx — 1) » 
f 3> \/ Probi,{hi), only a single term corresponding to 
/ii contributes significantly to the summation, and then 
the summation is nearly independent of t. It is clearly 
seen from Fig. 6, where Prob^Ti > t) is plotted t for 
L = 20 in the Id Oslo model. We can identify three steps 
in the plot. Each step in the curve can be associated with 
a unique value of hi (hi = 28, 29 and 30) and steps ap- 
pear at the corresponding value of Probi,(h\) along the 
y-axis. 

This explains the steps like structure of ProbL (2i > t) 
as a function of t. Also, the function decays roughly as 
1/t since we must have Probi,(h\) ~ 1/t for the term 
to contribute in Eq. 10. If h*(t) is the value of hi that 
contributes most in Eq. (JjJJ, the value of h\(t) is given 
by the condition Probt,(h\(t)) k, 1/t. Substituting this 
condition in Eq. (JTJ, we get 



Kit) 



L^/t 



(6) 



Thus the t dependence of h* (t) comes through the scal- 
ing variable tL^^ 1 = r (say). Then for r large the argu- 
ment establishing the 1/t dependence of Prob^iti > t) 
given above is quite robust. However more careful anal- 
ysis of Eq. (J5J shows that there is also a logarithmic 
multiplicative correction factor with the 1/t decay of 
Prob L (Ti >t). 

For large L and t, the terms, which contribute to 
ProbL(Ti > t) in Eq. JSJ, correspond to the values of 
hi for which hi <C hi- Substituting the scaling form 
of Prob L (h) = j^gi^jr^) (see Eq. Q) in Eq © and 
putting x = (h — hi)/L Ul , we get the following. 



Prob L (Ti > tL" 1 ) ~ / dxg(x)exp[-g(x)T} (7) 
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FIG. 7: The cumulative probability ProbL(Ti > i) for resi- 
dence time at the first site has been plotted against the scaled 
residence time t/L 0,25 for lattice sizes L = 300,400 and 500 
in the Id Oslo ricepile model. Total 10 7 grains were added. 



where t — t/L Ul . We have assumed that the probability 
distribution ProbLihi) <C 1, but is not rapidly decay- 
ing so that ProbLihi — 1) / 1 Probh(h\) ~ 1, and then the 
summation of Eq. © can be replaced by an integral over 
the scaled variable x. Actually in the real simulation (or 
experiment) for large L, this is the region of t we explore, 
we cannot go too far down the tail of ProbhiTi > t). 

Now it is easy to see the origin of logarithmic correction 
if we choose a particular form of the scaling function g(x) 
as exp(— |a:| a ) for i < 1, with a > 0, and try to find 
out the large t behaviour of the above equation in terms 
of the scaling variable r = t/L Ul . We first substitute 
s = exp(— \x\ a ) in Eq @ and get, 



Prob L {Ti > tL U1 ) 



ds 



-ln(s)] 



— exp(-sr) (8) 



The asymptotic behaviour of the of the above integral 
for r large is easy to evaluate, giving 



Prob L {Ti > tL u 



l/r(lnr) 



(9) 



In Fig. 7 we have plotted Prob^Ti > t) against scaled 
variable t/L Ul with loi = 0.25 for large values of system 
sizes L = 300, 400 and 500 in the Id Oslo model. We fit 
the scaled curves with a functional form given in Eq. l{9"|) 
with a = 2 since the scaling function g(x) is Gaussian 
near x = 0. 

Now it is clear that there is a logarithmic multiplicative 
factor in 1 /t decay and we take account of this logarith- 
mic multiplicative correction by writing the cumulative 
probability as given below. 



Prob L {Ti > t) ~ L Ul [tlog Sl (tL- U1 )] 



(10) 
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FIG. 8: The cumulative probability Probi,(Ti > t) for resi- 
dence time at the first site has been plotted against the scaled 
time t/L - 25 for lattice sizes L = 20,25,35 and 50 in the Id 
Oslo ricepile model. Total 10 9 grains were added. 



where we have used the fact that the answer is function 
of the scaling combination tL~ Ul and 8\ = (a — 1) /a. As 
a check we calculate the average residence time at site 1 
as given below. 



Hence the distribution of their sum T = ^ Tj would also 
be of same form. 

Even though the cumulative distribution of residence 
times Ti at any site i has steps like structure, the step- 
structure may be washed out in the sum J2i Ti- 

Results of the numerical simulation for distribution of 
the total residence times using a total 5 x 10 7 grains are 
shown in Fig. 9. We see that the steps are not seen in 
the distribution Pro6i(T > t) for different values of L 
for the range of the total residence times reached in the 
simulation (T < 10 8 ). The function Prob L (T > t) is 
much smoother than the function Probh{Ti > t). How- 
ever for small values of L (say for L < 20) , various curves 
of ProbL{T > t) still cross each other at large times. But 
for larger values of L, we don't see any inter-crossing of 
the curves in the times reached in our simulation (ex- 
cept at the tail where the data is less reliable due to the 
statistical fluctuations). 

In analogy with results for the the distribution of T\ , 
We can expect the behaviour of the cumulative distri- 
bution Prob^T > t) to be a scaling function of t/L u , 
where the exponent u> is different from uji defined earlier. 
So we write 



Prob L (T ><) = /(- 



(12) 



Ti,„ 



Prob L (Ti > t)dt ~ ]n(T l max ) 1 ' 



Si 



(11) 

The upper cutoff on the timescale is provided by 
1/ Prob^ih = h m i n ), which is the average time interval 
between successive returns to the minimum height. As- 
suming the scaling function g{x) varies as exp(— |a;| a ) for 
x <C — 1 and then putting T\^ max — exp[fcT Q ( 1_Wl )] in the 
above equation, we see that (Ti) is proportional to L. 

For the Id Oslo model, numerical estimate from the 
simulation gives lo\ w 0.25. Assuming the value a ks 
4 (argued in section 5), we get Si fts 0.75. In Fig. 8 
we have plotted ProbL^Ti > t) versus a scaled variable 
t/L" 1 where uj x = 0.25 for L = 20,25,35 and 50 and 
fitted the envelop formed by steps in the curves with a 
function l/x(lnx) Sl where Si = 0.75. We see that we get 
a reasonable fit to the data. We note that the multiplying 
logarithmic factor is necessary to get a good fit to 1/t 
dependence. 



where the scaling function f(x) varies as 1/ [a: (log a:)" 5 ] 
for large x, and the exponent S would also be different 
from Si defined earlier. Using the condition that the 
mean residence time in the pile is equal to the mean 
active mass in the pile, and hence scales as L 2 , can be 
used to determine S in terms of u> and 7 by integrating 
ProbhiT > t) over t upto the cut-off time scale exp(/cL 7 ). 
Now we get, 



*=l-(2-w)/7 



(13) 



In Fig. 10 we have plotted ProbL(T > t) versus scaled 
variable t/L u where u w 1.25 llO(. We get the value of 
S approximately equal to 0.81 from Eq. (T3J, assuming 
7 = 3 (argued in section 5). The fit is seen to be very 
good. In the numerical analysis of Christensen et. al. |5j, 
no logarithmic factor was used, and the data was fitted 
with a larger effective exponent, i.e. , 1/t 1 ' 22 decay. 



4. THE DISTRIBUTION OF RESIDENCE TIMES 
T IN THE OSLO MODEL. 

The arguments given in the previous section are easily 
extended to distribution of the residence times Ti with 
i =/= 1, and we conclude that they would also have a simi- 
lar 1 [t distribution with same logarithmic correction fac- 
tor which is for Ti , so long as i is not near the right end. 



5. PROBABILITY OF MINIMUM SLOPE IN THE 
OSLO MODEL 

The function Pro6i(/ii) can be exactly calculated nu- 
merically for small L using the operator algebra satisfied 
by addition operators 0|- We denote any stable config- 
uration by specifying slope values at all sites from i = 1 
to i = L, e.g., 1 122 21). Whenever slope zi becomes 2 
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FIG. 9: The cumulative probability ProbL(T > t) versus time 
t for lattice sizes L = 8, 11, 15 and 20 in the Id Oslo ricepile 
model. Total 5 x 10 7 grains were added. 
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FIG. 11: Probability of occurrence of minimum slope config- 
uration, calculated exactly, is plotted versus system sizes L in 
the Oslo ricepile model in the semi log scale. We calculated 
for L = 1 to 12. 
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FIG. 10: Scaling collapse of Prob^iT > t) versus scaled res- 
idence time t/L 1 ' 25 in the Id Oslo ricepile model for lattice 
sizes L — 20, 25, 35 and 50. Total 5 x 10 7 grains were added. 



after additions or toppling at site i, we denote such slope 
by 2, i.e., |...2...). "Bar" denotes that the site may topple 
or become stable with probability q or p respectively. 

Using these two toppling rules repeatedly and the 
abelian property of the Id Oslo ricepile model, we can 
relax any unstable configurations. For example, if we 
relax 1 22) for L = 2, we get the following sequence, 
|22) ->p|22) + g|12) ^ p 2 \22) +pq\l2) +pq\l2) + q 2 \2l) -> 
. . . -f p 2 \22) + (p + p 2 )q\12) + (p + p 2 )q 2 \21) + (p + 
f> 2 )g 3 |02) + (l+p)g 4 |ll). 

The probability of maximum slope configuration (i.e., 
when hi = 2L) can be easily calculated. We start with 
the unstable configuration |22...2). The probability that 
no site topples in this unstable configuration is p L and 
this is the probability of the maximum slope configu- 
ration (i.e. hi = 2L) in the steady state. That this 



probability varies as exponentially with L can be incor- 
porated in the scaling hypothesis by assuming that the 
scaling function g(x) in Eq. Q varies as exp(— aa; 1_ "i ) 
for x 3> 1 where a is a constant. 

The probability of the minimum slope configuration 
cannot be calculated so easily. However we argue below 
that this probability asymptotically varies as exp(—nL 3 ) 
where k is a constant. 

Firstly, the above calculation for L — 2 showed that, in 
the steady state, the probability of the minimum slope 
configuration is Q(q 4 ). For L = 3, we calculated this 
probability explicitly |l9| which is 0(q 10 ). Similar anal- 
ysis, for other values of L = 1 to 20, shows that the 
probability of minimum configuration is 0(q m ^), where 
m L is exactly given by the formula L(L+l)(L+2)/6. The 
coefficient of q m ^ in the probability is harder to compute 
explicitly for large L. We conjecture that this simple for- 
mula holds true for all L. Then for sufficiently small q, 
the probability of minimum height configuration in the 
lei Oslo model varies as exp[— n(q)L 3 ], where n(q) is a 
g-dependent function. Then, as there is no change in the 
behaviour of the Oslo ricepile expected, as a function of 
q, this behavior should persist for all non-zero q. For 
the scaling function, this would imply that g(x) varies as 

3 

exp[— k(<7)|x| ] for x 3> 1. 

We have calculated, ProbL(slope = 1), i.e. , the prob- 
ability of the minimum slope configuration, exactly nu- 
merically for q = 0.50, 0.60, 0.75 for L = 1 to 12 and the 
logarithm of Prob^slope — 1) has been plotted versus 

L(L + l)(L+2) • T7- 11 

*■ q ' m lig. 11. 

More specifically, consider a very low-slope unstable 
configuration 1 1 1 ... 12) which has total L number of grains 
with slopes 1 at all sites except at the last site with slope 2 
and estimate the probability to go to the minimum slope 
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FIG. 12: The residence time distribution ProbL{T l t > i) of 
grains at the corner site versus time t in model-B for lattice 
sizes L = 12, 15 and 20. Total 10 6 grains were added. 



FIG. 13: The cumulative probability distribution function 
Probh{T l 1 > t) versus scaled scaled residence time T l j/i ' 2 
in model-B for lattice sizes L = 12, 15 and 20. Total 10 6 
grains were added. 



from this configuration. To do this, we have to remove L 
grains, and each grain has to be moved a distance of 0(L) 
on the average. Thus we need 0(L 2 ) steps for large L, 
and each step requires a factor q in probability. Actually 
the probability of transition from this configuration with 
height h m i n + 1 to the minimum slope configuration (with 

LjL + l) L(L + 1) 

height hmin) is 0(q 2 ) and the coefficient of q 2 
in this case is exactly 1. Now the probability of minimum 
slope can be written in a general form as given below. 

Prob{slope — 1) ~ cxp[- n(q).L 3 ] (14) 

where k(0) = 00 and k(1) = 0. The different asymp- 
totic behaviour of large deviations in g{x) is somewhat 
unexpected, but has been seen in other problems, such 
as distribution of the large deviation of current in the 
asymmetric exclusion process in a ring |20l| . 
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FIG. 14: The residence time distribution Probh{T > t) in 
model-B for lattice sizes L = 12, 15 and 20. Total 10 6 grains 
were added. 



6. GENERALIZATION TO OTHER MODELS 

In this section we present the simulation results in 
other models and show that the cumulative distributions 
Probh(T\ > t) and Probi,{T > t) have same \ jt power 
law behaviour for large t, but with different logarithmic 
corrections. 

Model-B : Ricepile model in two dimension. 

Now we present the simulation results for 2d ricepile 
model. We add marked grains at the corner site, i.e., at 
(1,1). We simulated this model choosing p = 0.75 and 
q = 0.25, and study the residence time distribution of 
grains at the corner site (1,1). The standard deviation 



cr/j i of height h 1 1 at the corner site about the mean 
varies as L Wl where we estimated u>i ~ 0.2 from the sim- 
ulation. We added total 10 6 grains. 

We have plotted various curves for cumulative distri- 
bution function Prob^^Tx^ > t) of residence time T^i 
at the corner site versus times t for L = 12, 15 and 20 
in Fig. 12 and we see steps like structure appearing for 
t > 50. Various curves for different L inter-cross each 
other many times as seen in the Id Oslo ricepile model. 
In this case also, steps like structures are not log pe- 
riodic as the step-length in each curve increases on log 
scale when going down the curve. Any simple finite size 
scaling does not work as in the case of the Id Oslo ricepile 
model. However in the Fig. 13 we plotted various cumu- 
lative distributions ProbL,{T\ t \ > t) versus a scaled time 
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FIG. 15: Scaling collapse of various ProbhiT > t) versus 
scaled variable t/L 2 ' in model-B for lattice sizes L — 15,20 
and 27. Total 10 6 grains were added. 
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FIG. 16: The cumulative probability ProbL(T\ > i) versus 
time Ti for lattice sizes L = 25, 50 and 100 in model-C. Total 
10 6 grains were added. 



t/L" 1 with u>i = 0.2. We see that decay of the envelop 
formed by various steps in different curves fit well with 
the function 2/[x(lnx) ' 5 ] where the logarithmic correc- 
tion factor is according to Eq. ©. 

Using a similar argument to the Id Oslo model, in this 
case, we must have 7 = 4. To get to the minimum slope 
configuration, we will have to topple 0(L 3 ) grains and 
each grain O(L) times. As the average mass of the pile, 
in this case, varies as L 3 , Eq. 113(1 is modified as given 
below. 

6 = l-(3-u)/ r 

In Fig. 14 we have plotted various distribution of the 
total residence time, Probh(T > t), versus time t for dif- 
ferent L — 12, 15 and 20. In Fig. 15 we have plotted 
ProbL(T > t) for different L against the scaling vari- 
able t/L'* 1 where u) ps 2.0. Now we can estimate 6 to 
be approximately 0.75 from the above equation. In Fig. 
15 we fit the scaling function for Prob^(T > t) with 
0.2/[x(lnx) a75 ] which seems to be a reasonable fit. 

Model-C : The Local limited model. 

Residence time T\ : Since in model-C, grains are 
added randomly everywhere in the pile, average resi- 
dence time (Ti)/ n of a grain added at height h\ varies as 
l/[p\ProbL{hij\ according to Eq. JIJ|. Now there will be 
an extra l/L factor inside the exponential in the Eq. J7J . 
Consequently the scaling variable r = in Eq. (JJJ is 
replaced by r = L i* ui and Eq. 10 is modified to 

Prob L (Ti > tL x+ ^) ~ (far) ~ (15) 



Similarly average residence time at the first site equals 
to (hi)/pi (proof is similar as for (Ti)^ in Eq.^J which, 
in this case, varies as L 2 . This can be checked directly 
by integrating the above equation upto the cutoff time 
scale as done in Eq. i|ll|) . 

Total number of grains added in the pile are different 
for different L so that 10 5 grains are added at the first 
site. The standard deviation of height fluctuations 
at the first site varies as L" 1 with system size L, where 
oji w 1/3 UK. We have plotted Prob L (T x > t) for differ- 
ent values of L in the log-log scale in Fig. 16. We note 
that, unlike in the Id Oslo ricepile model, the cumulative 
probability here is smooth (except at the tail due to sta- 
tistical fluctuations) and monotonic function of L for a 
fixed t. This is due to the fact that the probability distri- 
bution ProbL(hi) of height at first site is not as sharply 
decaying function for hi hi as it was in the Id Oslo 
model. In fact, in Fig. 17 we get a good scaling collapse 
of various Prob^Ti > t) for different L using the scaled 
residence time t/L 1+UJl where uji ps 1/3. The scaling 
function is fitted well with the function l/[x(lnx)(" _1 ) // "] 
for x > 1, taking a = 2 (see Eq. ©). 

Total residence time T : In Fig. 18 we have plot- 
ted various Prob^i/T > t) versus residence time t for 
lattice sizes L — 50, 70 and 100. Total 10 6 grains were 
added in this case. In Fig. 19 we have collapsed vari- 
ous ProbhiT > t) for different L using the scaled vari- 
able t/L" where uj ps 1.5. We fit the scaling function of 
the cumulative distribution with l/2x[ln(2a;)]' 5 for x ^ 1 
where 6 ps 0.63. 
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FIG. 17: Scaling collapse of Probi,(Ti > i) for lattice sizes 
L = 25, 50 and 100 in model-C. Total 10 6 grains are added. 



FIG. 19: Scaling collapse of ProbijiT > t) versus scaled time 
t/L 1 ' 5 for lattice sizes L = 50, 70 and 100 in model-C. Total 
10 6 grains were added. 
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FIG. 18: The distribution function ProbL(T > t) versus time 
t in model-C for lattice sizes L = 50, 70 and 100. Total 10 6 
grains were added. 
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FIG. 20: The cumulative probability Pro6i(Ti > t) versus 
time t for lattice sizes L = 75, 100 and 130 in model-D. 10 4 
grains were added at the first site. 



Model-D : Variation of the Local Limited model. 

In the model-D, the standard deviation <jf ll of the 
hight fluctuation at site 1 scales with L as L" 1 where 
we found u>i f» 0.53. In Fig. 20 we have plotted various 
Prob^{T\ > t) against the residence time t at the first 
site for L = 75, 100 and 130. We added 10 4 grains at 
the first site. In Fig. 21 we have plotted Prob^lTx > t) 
versus scaled time t/L 1+UJl with w\ ~ 0.53 and get a 
good scaling collapse of all the curves for various L. We 
fit the scaling function with 0.73/[0.5a;(ln0.5x) (Q " 1)/Q ], 
using Eq. and putting a — 2, as done in the model-C. 

In Fig. 22 we have plotted various ProbL(T > t) for 
the total residence time versus t for lattice sizes L = 
50, 70 and 100. Total 10 5 grains were added in this case. 
In Fig. 23 we have plotted Pro6i(T > t) versus scaled 



time t/L" where u> ~ 1.5 and we get a good collapse for 
the scaling function which fits reasonably well with the 
function O^/^ln^)]' 5 for S w 0.63. 



11. SUMMARY AND CONCLUDING REMARKS 

To summarize, in this paper, we studied distribution 
of the residence times of grains in various ricepile mod- 
els. We reduced the problem of finding the residence time 
distribution of grains at a particular site to that of de- 
termining the distribution of first return time of height 
at the site to the same value. The result that the prob- 
ability of the residence times Ti at site i or the total 
residence time T in the pile, being greater than or equal 
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FIG. 21: Scaling collapse of Probi,(Ti > i) for lattice sizes 
L = 70, 100 and 130 in model-D. 



FIG. 23: Scaling collapse of ProbL(T > t) against the scaled 
variable t/L 1 ' 5 for lattice sizes L — 50, 70 and 100 in model-D. 
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FIG. 22: The cumulative probability Profit (T > t) versus 
time t for lattice sizes L = 50, 70 and 100 in model-D. Total 
10 J grains were added. 



size L is possible when the probability distribution func- 
tion ProbL(hi), where hi is the height at site i, sharply 
decays for hi < hi. However this non-monotonicity is 
seen only for t 3> t*(L) where t*(L) increases with in- 
creasing values of L, and hence may be harder to observe 
in real experiments. 

It is important to note that if we change the transfer 
rule of grains, the distribution of residence times may 
change completely. The rule, chosen in this paper is 
called first-in-last-out rule. We may employ some other 
rules, such as the first-in-first-out rule, adding the grain 
at the top of the stack but take out grains from the bot- 
tom of the stack. A different rule would be to add and 
take out grain from a stack in random order. In these 
cases, there are no sites with deeply buried grains and 
the residence time distribution will be similar to that in 
the critical height models studied earlier by us. 



to t, decay as power law 1/t is valid for a large class of 
sandpile models, where height fluctuation at a particular 
site grows with the system sizes, and is independent of 
dimensions. It depends only on the fact that there are 
some deeply buried grains which come out only in rare 
fluctuations, i.e. , slope of the pile becomes very close to 
the minimum slope. It is important to note that, since 
the total residence time T is sum of TVs, the probability 
of T = is very small, and clearly our analysis can- 
not predict the behaviour of the cumulative probability 
Prob L (T > t) for t < L u . 

We also found that cumulative probability Probi,{T\ > 
t) is non-monotonic with system size L for any fixed 
t for some of the ricepile models. The non-monotonic 
behaviour of the cumulative probability distribution 
Probj^^Ti > t) of residence times at site i with system 
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